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ABSTRACT 

We study the impact of inhomogeneous hydrogen reionization on the thermal evolution of the 
intergalactic medium (IGM) using hydrodynamic + radiative transfer simulations where reionization 
is completed either early (z ~ 9) or late (z ~ 6). In general, we find that low-density gas near 
large-scale overdensities is ionized and heated earlier than gas in the large-scale, underdense voids. 
Furthermore, at a later time the IGM temperature is inversely related to the reionization redshift 
because gas that is heated earlier has more time to cool through adiabatic expansion and Compton 
scattering. Thus, at the end of reionization the median temperature-density relation is an inverted 
power-law with slope 7— 1 ~ —0.2, in both models. However, at fixed density, there is up to order unity 
scatter in the temperature due to the distribution of reionization redshifts. Because of the complex 
equation-of-state, the evolved IGM temperature-density relations for the redshift range 4 < z < 6 
can still have significant curvature and scatter. These features must be taken into account when 
interpreting the Lya absorption in high redshift quasar spectra. 

Subject headings: cosmology: theory - large-scale structure of universe - intergalactic medium - 
methods: numerical - hydrodynamics - radiative transfer 



1. INTRODUCTION 

The cosmic reionization of hydrogen fundamentally 
changes the thermal and ionization conditions for the 
high redshift (z > 6) intergalactic medium (IGM), con- 
verting a cold and neutra l gas into a warm and highly ion- 
ized one (see |Loeb|2008| for a review) . During this inho- 
mogeneous process, the radiative transfer (RT) of the UV 
field, primarily generated by stellar sources, will proceed 
such that large-scale, overdense regions near sources are 
generally ionized and heated earlier th an large-scale, un- 



derdense regions far from sources (e.g. Barkana fc Loeb 



2004 
et al. 



Wyithe fc Loe b 2004; Furlanetto et al.||2004|]Lee| 
2008 f. The photo-ionization of hydrogen (HI — ► 



Hll) and helium (Hel — > Hell) is expected to heat the 
IGM to temperatures of > 10 4 K, with higher tempera- 
tures obtained for harder radiation spectra (e.g.|Miralda-| 
Escude fc Ostriker||1990] |Miralda-Escude fc Rees||1994p . 
Thereafter, the thermal evolution will be governed by 
adiabatic cooling/heating for underdense/overdense gas, 
Compton scattering on the Cosmic Microwave Back- 
ground (CMB), atomic line cooling, shock-heating, and 
additional photo-heating only where recombination is ef- 
ficient. 

In principle, the temperature of the IGM can be quan- 
tified through the thermal broadening of the Lyman al- 
pha (Lya) forest lines (e.g. [Schaye et al. | |1999[ |2000| [RT 
cotti et al.||2000| |McDonald et al.||200ip tmd feeJeans 
smoothing of the L ya forest fl ux power spectrum (e.g. 



Theuns et al. 2000 Zaldarriaga et al. 2001 Viel et al 



2004p . Measurements of the temperature can then be 
used to put constraints on when t he reionization of try 



drogen and helium o ccurred (e.g. Theuns et al. 2002 



Hui & Haiman||2003l . Note that the full reionization o 
I] 



helium (Hell 
shifts (z ~ 3) 



Helll) is expected to occur at lower red- 



(see 


Furlanetto & Oh 2008 


Bolton et al.|2008 


McQuinn 


et ai 


2008 


for recent modelling). The current observa- 
;raints have large uncertainties, although bet- 


tiona 


1 cons 



will be required to interpret them. 

In this Letter we focus on the impact of inhomogeneous 
hydrogen reionization on the thermal evolution of the 
IGM. In particular, we study the equation-of-state and 
measure the temperature-density relation and its scatter 
using hydrodynamic + RT simulations. |Hui fc Gnedin 
( 1997 1 have studied the effects of uniform reionization 
and found that the temperature-density relation is well 
approximated by a positive power-law slope for the low- 
density IGM. However, inhomogeneous reionization can 
give rise to a much different equation-of-state, especially 
right after completion. Specifically, the large-scale, un- 
derdense voids are initially hotter than higher density gas 
near sources, that having been ionized and heated earlier 
have also had more time to cool. This gives rise to an 
equation-of-state characterized by an inverted (negative 
slope) temperature-density relation. At later times, the 
underdense gas cools adiabatically faster and the tem- 
perature contrast decreases and eventually reverses back 
to the positive slope. 

2. SIMULATIONS 

We adopt a hybrid approach in modelling cosmic reion- 
ization that allows us to cover a larger dynamic range and 
satisfy the requirements of having sufficiently high reso- 
lution to capture small-scale structure, and a sufficiently 
large volume to reduce sample variance ( Barkana & Loeb 
20041. Below, we summarize two recent numerical sim- 
ulations incorporating N-body, hydrodynamic, and RT 
algorithms to solve the couple d evolution of the dark 
matter, baryons, a nd radiation ( |Trac fc Pen||2004 2006 
Trac & Cen 20071. The simulations are based on the 



ter measurements will come in time and improved models 



recommended WIVIAP 5- year cosmological parameters: 
Cl m = 0.28, n A = 0.72, n b = 046, h = 0.70, a 8 = 0.82, 
and n s = 0.96 flDunkley etltil2008| ) . 

Our hybrid approach is divided into two major compo- 
nents. The first task involved running a high-resolution 
N-body simulation with 3072 3 dark matter particles on 
an effective mesh with 11520 3 cells in a comoving box, 
100 Mpc hr 1 on a side. We identified collapsed dark mat- 
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Fig. 1. — Low-resolution images of a slice of size (lOOMpch -1 ) 2 with a thickness of two hydro cells (130kpc/i — 1 ) from the late reionization 
model. Left: Redshift at which Jeans-size elements of the IGM get reionized. Right: Temperature field at the end of reionization. In general, 
low-density gas that is ionized and heated earlier has lower temperature at late times. 

ter halos on the fly using a friends-of-friends algorithm, 
with a linking length b = 0.2 times the mean interparticle 
spacing, in order to model radiation sources and sinks. 
With a particle mass resolution of 2.68 x 10 6 M ft -1 , 
we can reliably locate all dark matter halos with virial 
temperatures above the atomic cooling limit (T ~ 10 4 
K) with a minimum of 



40 particles (Hcitmann ct al. 
20071, and half of this collapsed mass budget is resolved 
with > 400 particles per halo. Our halo mass functions 
are in very good agreement with other recent work (e.g 



Reed et ai.|2007l|Lukic et al |2007| |C~ohn fc White|2008 |. 

The second component of the work produces a series 
of hydrodynamic + RT simulations with moderate reso- 
lution, but incorporating subgrid physics modelled using 
the high-resolution information from the large N-body 
simulation. For the first two simulations in the series, we 
consider basic models where reionization is completed 
early (z ~ 9) and late (z ~ 6) to roughly reflect our cur- 
rent state of knowledge on the details of the reionization 



process (e.g. |Cen||2003| |Wyithe fc Loeb|2003| |Haiman fc 



Holder 2003| |Fan et al.||2006ap . Radiation sources are 
prescribed and star for mation rates calcul ated using the 
halo model described in |Trac fc Cen (2007). Her e we con- 
sider only Population II stars from starbursts ( Schaerer 
2003 1 as contributing to the ionizing photon budget. We 



also neglect subgrid clumping and self-shielding of dense 
absorbers such as minihalos or Lyman limit systems, but 
will study their effects in ongoing work. 

Each simulation utilizes equal numbers (N — 1536 3 ) of 
dark matter particles, gas cells, and adaptive rays, where 
for the latter, we track 5 frequencies above the hydrogen 
ionizing threshold of 13.6 eV. The photo-ionization and 
photo-heating rates for each cell are calculated from the 
incident radiation flux and used in the non-equilibrium 
solvers for the ionization and energy equations. The 
initial conditions are generated with a common white 
noise fiel d and a linear tran sfer function calculated with 
CAMB (Lewis et al. 20001. The simulations were run 



using the NASA Columbia Supercomputer. 

3. RESULTS 

3.1. The Redshift of Reionization 

Over the course of each simulation we track when a 
gas cell first becomes more than 99% ionized, enabling 
us to construct a 3D reionization-redshift field z rc i on (x) 
in parallel to the gas density p(x) and temperature T(x) 
fields. In Fig. [I] we show a slice from the late reioniza- 
tion simulation box to illustrate the inhomogeneity of the 
reionization process and the large-scale correlations. The 
fields have been smoothed on cells of comoving length 
130 kpc fe -1 , which is close to the Je an's or filtering scale 
( |Gnedin fc Hui||1998| lGnedin||2000| ). In general, z roion is 
highly correlated with the large-scale density field, as 
HII regions originate around biased, overdense sources 
and final overlap occurring in the large-scale, undcrdense 
voids. We find that the strong and positive correlation 
extends down to scales ~ 1 Mpc hr 1 . 

The temperature field at the end of reionization shows 
considerable complexity with some interesting character- 
istics. First, we note the presence of cool (T < 10 4 K), 
low-density gas just outside the shock-heated filamentary 
and halo gas. Second, the low-density gas gets progres- 
sively warmer (T > 10 4 K) as one moves away from 
sources and towards the voids. It is apparent from the 
maps that the temperature is inversely related to z ro i on 
for the low-density IGM and we will quantify this relation 
in the next section. 

3.2. The Temperature- Density Relation 

Fig. [2] shows the gas temperature-density relations at 
the endof reionization and at redshift z = 4 for the early 
and late reionization models. For each density bin of 
width 0.05 dex, we plot the median temperature and the 
2-<7 spread around it. At the end of reionization, the two 
models have complex but similar T-p relations, driven by 
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Fig. 2. — Gas temperature-density relations at the end of reionization (left) and at a redshift z = 4 (right) for the early (red, thick lines) 
and late (blue, thin lines) reionization models. For each density bin of width 0.05 dex, we plot the median temperature (solid lines) and 
the 2-cr spread (dashed lines). At the end of reionization, the two models have similar inverted power-law T-p relations for the low-density 
IGM, with slopes 7 — 1 ~ —0.2 and order unity scatter. By 2 = 4, the early model has settled into a rising power-law with only small 
scatter. However, the late model is closer to yielding an isothermal equation-of-state with still significant scatter of thousands of kelvins. 



photo-heating, shock-heating, and cooling. For the low- 
density IGM where photo-heating dominates over shock- 
heating, we find that the median temperature approxi- 
mately follows a power-law relation, 



T = T (1 + S) 



7-1 



(1) 



with a slope 7 — 1 ~ —0.2. The median temperature 
has an inverted form, but the difference in the median 
temperatures at 6 = and S = —0.8 is small compared to 
the scatter in the temperature for underdense gas. The 
order unity scatter in the temperature at a fixed density 
is explained below. 
We first note that the IGM density on the Jean's scale 
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Fig. 3. — End-of-reionization temperature as a function of the 
reionization redshift for low-density (<5 = —0.5 ± 0.05) gas, for the 
early (red, thick lines) and late (blue, thin lines) models. The 
two contours specify the 1-u (solid lines) and 2-<r (dashed lines) 
boundaries for the cumulative probability distributions. 



is only weakly correlated with the density averaged on 
large scales, which are still in the linear regime, low- 
density gas is also found close to large-scale overdense 
regions around sources and it preferentially gets ion- 
ized and heated earlier than gas at the same density 
in the large-scale voids far away from sources. Once 
heated, the gas cools through adiabatic expansion and 
Compton scattering. There is only minor additional 
heat input from subsequent photo-ionization of the small 
residual neutral fraction because the recombination time 
is longer than the Hubble time at low densities. In 
Fig. [3] we plot the temperature of gas cells with den- 
sities S = — 0.5 ± 0.05, taken at the end of reionization, 
as a function of z rc ion- Indeed, we find a strong correla- 
tion whereby low-density gas with earlier z rc ion has lower 
temperatures. Thus, the variations in the temperature 
at a fixed density in Fig. [2] are simply correlated with 
the timing of reionization. We expect a more negative 
slope 7 — 1 and a larger spread in temperature for more 
a extended reionization epoch. 

Fig.[2]also shows the evolved temperature-density rela- 
tions at z = 4. For the early model, the scale factor has 
roughly doubled since the end of reionization and enough 
time has passed for the photo-ionized IGM to cool and 
settle into a monotonically rising curve. The T-p relation 
for underdense gas can be approximated by a power-law 
with an effective slope 7 — 1 ~ 0.5, which is close to the 



corresponding value calculated by Hui & Gnedin (1997 
see Fig. 4) for uniform reionization and smaller than the 
maximum value of 0.62 obtained when far after reion- 
ization. The adiabatic and Compton cooling processes 
have also significantly reduced the scatter such that the 
96% spread is now less than a thousand degrees K, also 
similar to what is found for uniform reionization. 

However, the T-p relation for the late model is still 
relatively complex, with the lowest density gas having 
an isothermal equation-of-state while near mean density, 
the effective slope is 7 — 1 ~ 0.47. Furthermore, the scat- 
ter is still significant with a 96% spread of several thou- 
sand degrees K. These characteristics are remnants of the 
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complicated conditions at the end of reionization. For 
uniform reionization, the scatter is several times smaller, 
approximately a thousand degrees K, and the curvature 
is also smaller, but still persistent since not enough time 
has passed for cooling to set up the power-law relation 



(e.g. Hui fc Gnedin 1997 



rather than the curvature 



Trac||2004| ). It is the scatter 
where we find the dominant 
difference because at the end of reionization, the inver- 
sion of the T-p relation is weak and not far from isother- 
mal. 

Interestingly, we find that the T-p relation for the early 
model at z = 6 (not shown in figures) is similar to that of 
the late model at z = 4. Thus, even for early reionization 
models that end at z ~ 9, the inhomogcncity of hydrogen 
reionization may still be imprinted in the spectra of high 
redshift quasars. Since the scatter in temperature is not 
random, but correlated with the redshift of reionization, 
which in turn is correlated with the large-scale density 
field, it may be possible to probe the correlated inhomo- 
geneity with the Lya fores t through statistics su ch as the 
flux power spectrum (e.g. Fang & White 20041. 



both the slope and scatter will depend on the reioniza- 
tion history, especially its duration. Furthermore, it is 
known that radiative transfer effects from Lyman limit 
systems can modify th e temperature distribution (e.g. 
Abel fc Haehnelt|1999[ ), although more for helium rather 
than hydrogen reionization because for the latter, the UV 
spectrum from stellar sources is relatively soft and only 
weak spectral filtering can occur. 

We conclude that at the high redshift range 4 < z < 6, 
the equation-of-state of the IGM is more complicated 
than the commonly assumed forms (isothermal or a tight 
power-law relation). It is important to keep this result 
in mind whe n interpreting the Lya absorption in quasar 



spectra (e.g. |Fan et al. 2002 2006b Lidz et al.||2006 



4. CONCLUSIONS 

We have studied the photo-ionization and photo- 
heating of the IGM from hydrogen reionization using 
two hydrodynamic + RT simulations. We considered 
two basic models in which reionization is completed early 
(z ~ 9) and late (z ~ 6) and found for both models that 
the temperature of a low-density region at a later time 
is inversely related to the redshift of reionization of that 
region, low-density gas near large-scale overdensities is 
ionized and heated earlier than gas in the large-scale, 
underdense voids. As a result, at the end of reionization 
the median temperature-density relation is an inverted 
power-law with slope 7 — 1 ~ —0.2. There is up to order 
unity scatter in the temperature at fixed density due to 
the distribution of reionization redshifts. We expect that 



Gallerani et al.|2008[ ). In fact , a recent analysis by Becker 
et al.| ( |2007| has already suggested that adopting an in- 
verted temperature-density relation instead of an isother- 
mal model may have profound implications on the inter- 
pretation of the reionization process based on the high 
redshift Lya forest. We will study the observational sig- 
natures of inhomogeneous hydrogen reionization in an 
upcoming paper. 



We thank C.-A. Faucher-Giguere, L. Hernquist, 
A. Lidz, M. McQuinn, and M. Zaldarriaga for many stim- 
ulating discussions. We also thank J. Chang at NASA 
for invaluable supercomputing support, N. Gnedin for 
his compilation of ionization and recombination rates, 
and D. Schaerer for the Pop II SEDs. This work is 
supported in part by NASA grants NNX08AL43G and 
NNG06GI09G, NSF grant AST-0407176, FQXi, and 
Harvard University funds. Computing resources were in 
part provided by the NASA High-End Computing (HEC) 
Program through the NASA Advanced Supercomputing 
(NAS) Division at Ames Research Center. 



REFERENCES 



Abel, T., & Haehnelt, M. G. 1999, ApJ, 520, L13 
Barkana, R., & Loeb, A. 2004, ApJ, 609, 474 

Becker, G. D., Rauch, M., & Sargent, W. L. W. 2007, ApJ, 662, 72 
Bolton, J. S., Oh, S. P., & Furlanetto, S. R. 2008, ArXiv e-prints, 
807 

Cen, R. 2003, ApJ, 591, 12 

Cohn, J. D., & White, M. 2008, MNRAS, 385, 2025 

Dunkley, J., Komatsu, E., Nolta, M. R., Spergcl, D. N., Larson, 

D., Hinshaw, G., Page, L., Bennett, C. L., Gold, B., Jarosik, N., 

Weiland, J. L., Halpern, M., Hill, R. S., Kogut, A., Limon, M., 

Meyer, S. S., Tucker, G. S., Wollack, E., & Wright, E. L. 2008, 

ArXiv e-prints, 803 
Fan, X., Carilli, C. L., & Keating, B. 2006a, ARA&A, 44, 415 
Fan, X., Narayanan, V. K., Strauss, M. A., White, R. L., Becker, 

R. H., Pentericci, L., & Rix, H.-W. 2002, AJ, 123, 1247 
Fan, X., Strauss, M. A., Becker, R. H., White, R. L., Gunn, J. E., 

Knapp, G. R., Richards, G. T., Schneider, D. P., Brinkmann, J., 

& Fukugita, M. 2006b, AJ, 132, 117 
Fang, T., & White, M. 2004, ApJ, 606, L9 
Furlanetto, S. R., & Oh, S. P. 2008, ApJ, 681, 1 
Furlanetto, S. R., Zaldarriaga, M., & Hernquist, L. 2004, ApJ, 613, 

1 

Gallerani, S., Ferrara, A., Fan, X., & Choudhury, T. R. 2008, 

MNRAS, 386, 359 
Gnedin, N. Y. 2000, ApJ, 542, 535 
Gnedin, N. Y., & Hui, L. 1998, MNRAS, 296, 44 
Haiman, Z., & Holder, G. P. 2003, ApJ, 595, 1 
Heitmann, K., Lukic, Z., Fasel, P., Habib, S., Warren, M. S., White, 

M., Ahrens, J., Ankeny, L., Armstrong, R., O'Shea, B., Ricker, 

P. M., Springel, V., Stadel, J., & Trac, H. 2007, ArXiv e-prints, 

706 

Hui, L., & Gnedin, N. Y. 1997, MNRAS, 292, 27 
Hui, L., & Haiman, Z. 2003, ApJ, 596, 9 

Lee, K.-G., Cen, R., Gott, J. R. I., & Trac, H. 2008, ApJ, 675, 8 



Lewis, A., Challinor, A., & Lasenby, A. 2000, ApJ, 538, 473 
Lidz, A., Oh, S. P., & Furlanetto, S. R. 2006, ApJ, 639, L47 
Loeb, A. 2008, ArXiv e-prints, 804 

Lukic, Z., Heitmann, K., Habib, S., Bashinsky, S., & Ricker, P. M. 

2007, ApJ, 671, 1160 
McDonald, P., Miralda-Escude, J., Rauch, M., Sargent, W. L. W., 

Barlow, T. A., & Cen, R. 2001, ApJ, 562, 52 
McQuinn, M., Lidz, A., Zaldarriaga, M., Hernquist, L., Hopkins, 

P. F., Dutta, S., & Faucher-Giguere, C. . 2008, ArXiv e-prints, 

807 

Miralda-Escude, J., & Ostriker, J. P. 1990, ApJ, 350, 1 
Miralda-Escude, J., & Rees, M. J. 1994, MNRAS, 266, 343 
Reed, D. S., Bower, R., Frenk, C. S., Jenkins, A., & Theuns, T. 

2007, MNRAS, 374, 2 
Ricotti, M., Gnedin, N. Y., & Shull, J. M. 2000, ApJ, 534, 41 
Schaerer, D. 2003, A&A, 397, 527 

Schaye, J., Theuns, T., Leonard, A., & Efstathiou, G. 1999, 

MNRAS, 310, 57 
Schaye, J., Theuns, T., Rauch, M., Efstathiou, G., & Sargent, 

W. L. W. 2000, MNRAS, 318, 817 
Theuns, T., Schaye, J., & Haehnelt, M. G. 2000, MNRAS, 315, 600 
Theuns, T., Schaye, J., Zaroubi, S., Kim, T.-S., Tzanavaris, P., & 

Carswell, B. 2002, ApJ, 567, L103 
Trac, H. 2004, PhD thesis, University of Toronto (Canada), Canada 
Trac, H., & Cen, R. 2007, ApJ, 671, 1 
Trac, H., & Pen, U.-L. 2004, New Astronomy, 9, 443 
— . 2006, New Astronomy, 11, 273 

Viel, M., Haehnelt, M. G., & Springel, V. 2004, MNRAS, 354, 684 
Wyithe, J. S. B., & Loeb, A. 2003, ApJ, 586, 693 
— . 2004, Nature, 432, 194 

Zaldarriaga, M., Hui, L., & Tegmark, M. 2001, ApJ, 557, 519 



